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Abstract 

In a previous paper the real evolution of the system of ODEs 

N 

1 m—1. m^n 

! z n =z n (t), z n = , n = l,...,N 

is discussed in Cat, namely the N dependent variables z n , as well as the N(N — 1) 
(arbitrary!) "coupling constants" g n m, are considered to be complex numbers, while 
the independent variable t ("time") is real. In that context it was proven that there 
exists, in the phase space of the initial data z n (0), i„(0), an open domain having 
infinite measure, such that all trajectories emerging from it arc completely periodic 
with period 2ir, z n (t + 2tt) — z n (t). In this paper we investigate, both by analyti- 
\Q | cal techniques and via the display of numerical simulations, the remaining solutions, 

and in particular we show that there exist many — emerging out of sets of initial 
data having nonvanishing measures in the phase space of such data — that are also 
completely periodic but with periods which are integer multiples of 2ir. We also elu- 
cidate the mechanism that yields nonperiodic solutions, including those characterized 
by a "chaotic" behavior, namely those associated, in the context of the initial-value 
problem, with a sensitive dependence on the initial data. 



> 

1 Introduction 

In a previous paper [1] one of us (FC) analyzed, in the complex domain, the dynamical 
system characterized by the Newtonian equations of motion 

JV 

Zn A~ Z n — ^ 9nm(z n 2 ffl ) , 

m=l, m^n 

z n = z n (t), z n = , n = l,...,N (1.1) 
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which obtain in the standard manner from the Hamiltonian 

1 N 1 N 

H(z,p) = -J2(pI + z 2 n ) + - 9nm(z n -z m ) 2 , (1.2a) 

rt=l n,m=l, n^m 

provided (as we hereafter assume — even though the main result of [1] holds without this 
restriction) 

9nm — 9mn- (■^•^v 

Here and below N is an arbitrary positive integer (N > 2), the indices n, m run from 1 to N 
unless otherwise indicated, underlined quantities are iV-vectors, say z = (z±,..., zjv), and 
all quantities (namely, the N "canonical coordinates" z n , the N "canonical momenta" p n , 
the N(N — l)/2 "coupling constants" g n m) are complex, while the independent variable t 
("time") is instead real. Of course the N complex equations of motion can be 

reformulated [1] as 2N real — and as well Hamiltonian [2] - - equations of motion, by 
introducing the real and imaginary parts of the coordinates z n , z n = x n + iy n , or their 
amplitudes and phases, z n = p n exp(i6 n ): 

N 

X n ~\~ X n = ^ ] r nm \_ a nm x nm { x n m ~ %nm) — ^nmUnm {llnm ~ ^nm)] 
m=l, mj^n 
N 

= r nm 1 9nm | COS (7 nm - 36» nm ), (1.3a) 

m=l, mj^n 
N 

jjn + Dn = Y] r nm [ a nmUnm (llnm ~ ^ x nm) ~ ^nmVnm { x nm ~ ^Unm)] 
m=l, m^n 
N 

= r nmbnm| sin(7„ m - 36»„ m ), (1.3b) 

m=l, mj^n 

where of course 

Zn — x n ~i~ ^Vni 

Znm = Z n Z m — X n X m + i{]Jn Urn) = x nm ~\~ iynm — ^nm 6Xp(i$ nm ), 

9nm — 0"nm ~i~ ib nm — |5»im| 6Xp(i7nm)' (1.3c) 

We shall return to the motivations for this choice to investigate the system 1)1. 1JI in the 
complex domain at the end of this Section 1 . 

If all the coupling constants coincide, g nm = g, the Hamiltonian system Q1.1JI is a well- 
known completely integrable many-body model (see for instance [2] and the references 
quoted there), and all its nonsingular solutions are completely periodic with period 2ir, or 
possibly an integer multiple of 2ir. (Indeed, in this integrable case the N coordinates z n (t) 
can be identified with the N zeros of a polynomial of degree N the coefficients of which are 
periodic in t with period 2ir, so that the set of these N zeros is also periodic with period 2tt, 
and each individual zero is therefore also periodic, although possibly with a larger period 
which is an integer multiple of 2ir due to a possible reshuffling of the zeros as the motion 
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unfolds; in the real case with all coupling constants equal and positive, g nm = g > 0, when 
the motions are confined to the real axis and no such reshuffling can occur due to the 
singular and repulsive character of the two-body forces, all real solutions are nonsingular 
and completely periodic with period 2ir, z(t+2n) = z(t); see for instance [2]). Here we focus 
instead on the more general case with completely arbitrary coupling constants g nm , which 
is generally believed not to be integrable. But even in this case — as proven in [1] — there 
does exist a domain of initial data 2.(0), i(0) having infinite measure in phase space - 
indeed, having a measure which is a finite fraction of that of the entire phase space - 
such that all the trajectories originating from it are completely periodic with period 2ir, 
z(t + 2ir) = z(t). As pointed out in [1], this is a somewhat surprising finding, inasmuch as 
it negates the expectation that, for a nonlinear dynamical system with several degrees of 
freedom that possesses completely periodic trajectories emerging from some specific initial 
data, any generic variation of these initial data destroy the complete periodicity of the 
trajectories or at least change their period — unless the system is essentially equivalent 
(say, via an appropriate change of variables) to a linear system (such as (jl.ip with all 
coupling constants vanishing, g nm = 0), which is certainly not the case for the system 
with arbitrary coupling constants g nm , at least not in any manner explicitly computable 
in closed form. 

But, as shown in [1], this fact is a rather elementary consequence of an approach 
(a "trick") introduced and rather extensively used recently to evince analogous results 
(see [2] and the references quoted there, as well as [3, 4, 5]). This same trick can as well 
be exploited to investigate the remaining solutions, namely those not belonging to the 
class of completely periodic solutions with period 2ir the existence of which was proven 
in [1]. This we do in the present paper, and we also confirm the insight thereby gained by 
exhibiting numerical solutions of performed via a computer code created by one of 
us (MS) [6]. In particular we demonstrate below the existence of open domains of initial 
data, having nonvanishing measures in the phase space of such data, which also yield 
completely periodic solutions, but with periods which are integer multiples of 2ir, and we 
also elucidate the mechanism that originates non-periodic and chaotic solutions. 

In the next Section (the first part of which is closely patterned after Section 2 of [1] and 
is reported here to make this paper self-contained) the "trick" mentioned above — which 
in fact amounts to a change of (dependent and independent) variables — is introduced. In 
Section 3 the analyticity properties in a time-like complex variable are investigated of the 
solutions of the system of ODEs, see (|2.5|) . related to via the trick. The implications 
of these findings as regards the periodicity of certain solutions of 1)1. 1|) are discussed in 
Section 4; in particular the mechanism that underlies the existence of completely periodic 
motions with higher periods (integer multiples of 2tt), and of non-periodic and chaotic 
motions, is elucidated. In Section 5 numerical examples of these trajectories are exhibited. 
In Section 6 some final remarks are proffered. Two appendices complete our presentation: 
Appendix A contains some developments (confined there not to interrupt the flow of 
presentation in Section 3) concerning the analytic structure of the solutions of the evolution 
equations (|2.5|) : Appendix B focuses on the three-body case (N = 3), in particular it reports 
its reducibility to quadratures and a discussion of the information about the analytic 
behavior of the solutions of (|2.5jl entailed by this fact. 

Let us end this introductory Section 1 with some remarks on the choice to investigate 
the motion determined by the Newtonian equations in the complex, rather than the 
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real, domain. A clear hint that, at least from a mathematical point of view, this is a more 
natural environment to work in, comes already from the treatment of in the inte- 
grable case with equal coupling constants, g nm = g, since, as mentioned above, it is then 
appropriate to identify the N particle positions z n = z n (t) with the N zeros of a time 
dependent (monic) polynomial of degree N in z, say Pn(z,z) such that Pw[t, z n (t)] = 
(see for instance Ref. [2]); and clearly the natural environment to investigate the zeros 
of a polynomial is the complex plane rather than the real line. In our context, an essen- 
tial motivation to work in the complex comes from the important role that analyticity 
properties play in our treatment, see below. Moreover motions roaming over the complex 
plane display a much richer dynamics than those restricted to the real line, especially in 
the case with singular interparticle forces, because of the possibility in the former case, 
but not in the latter, that particles go around each other. And it is then natural to 
re-interpret the (complex) iV-body problem as describing the (real) motion of N 

particles (in the plane), by introducing a one-to-one correspondence among the complex 
coordinates z n = x n + iy n , see (|1.3cJ) . and the real two- vectors in the plane r n = (x n ,y n ). 
But this approach, that is quite convenient to identify interesting many-body problems 
in the plane (see Chapter 4 of Ref. [2]), suffers in the present case from a drawback: the 
resulting many-body problem in the plane is not rotation-invariant, see (1.3). Indeed 
the many-body problem (1.3) is characterized by a (clearly rotation-invariant) harmonic- 
oscillator one-body force attracting every particle towards the origin, and by a (clearly 
not rotation-invariant) singular two-body force acting among each particle pair. As clearly 
seen from (1.3), the strength of the two-body force is proportional to the inverse cube of 
the (Euclidian) interparticle distance, hence it generally diverges when two particles col- 
lide; but this force also depends, both in modulus and direction, from the orientation of 
the interparticle vector, as well as from the phase 7„ m of the relevant interparticle coupling 
constant (jl.3c|) . Indeed the two-body force f nm acting on the n-th particle due to the ra- 
th particle is the two-vector f nm = r~^\g nm \(cos(-y nm - 39 nm ), sin(7„ m - 36 nm )), which is 
generally not aligned to the interparticle distance f nm = r nrn (cos9 nm ,sin9 nm ), see (1.3). 
The diligent reader is advised to try and become familiar with the specific implications of 
this fact, as they will be eventually helpful to understand the trajectories of the solutions 
of the model 1)1. lj) . see Section 5. 



In this section we describe the "trick" [1-4] that underlies our subsequent findings, and 
at the end we mention a remarkable property of the system of ODEs 1)1.1)1 . But firstly we 
rewrite, mainly for notational convenience, these equations of motion, 1)1.1)1 . as follows: 



2 The trick 



N 




(2.1a) 



m=l, mj^n 



and we note that the corresponding Hamiltonian reads 




(2.1b) 



71=1 



n,m=l, n^m 
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Here we introduce the additional constant u, which is hereafter assumed to be positive, 
lo > 0, and to which we associate the basic period 

T = 2tt/lu. (2.2) 

In the following it will sometimes be convenient to set u = 1 so that (|2.1aj) coincide with 
(jl.lj) and the basic period becomes T = 2ir, or to set instead u = 2ir so that the basic 
period becomes unity, T = 1. Of course these cases are all related via a rescaling of the 
dependent and independent variables, since clearly by setting 

z(t) = az(t), t = bt, u = u/b, g nm = [a 2 /b) 2 g nm (2.3a) 

with a, b two positive rescaling constants which can be chosen at our convenience, the 
ODEs (|2.1a|) get reformulated in a completely analogous "tilded" version, 

N 

z'n+ti 2 Z n = 9nm{z n — Z m ) 3 , (2.3b) 

771=1, mj^n 

where of course here the primes indicate differentiations with respect to the argument 
of the function they are appended to, z! = dz(t)/dt. Note in particular that by setting 
b = a 2 = uj one gets u = 1, g nm = g nm , namely the tilded version (j2.3b|) reproduces 
essentially (jl.lj) . 

Now, the "trick". Let us set 

z n (t) = exp(-iutX n (T), (2.4a) 
r = r(i) = [exp(2iwt) - l]/(2iw). (2.4b) 

As can be readily verified, this change of (dependent and independent) variables, (2.4), 
transforms (j2.1aj) into 

N 

c= E wcn-ur 3 . (2.5) 

m=l, mj^n 

Here and below appended primes denote derivatives with respect to the new independent 



variable r, while of course the dots in (j2.1aj) and below denote as usual derivatives with 
respect to the real time t. 

The change of variables (2.4) entails the following relations among the initial data, z(0), 
i(0), respectively C(0), C_{0), for (|2~Taj) respectively (|2~5|) : 

Zn(0) = Cn(0), (2.6a) 
z n (0) = C(0) - ^Cn(O). (2.6b) 

We now observe that, as the (real, "physical time") variable t varies from to T/2 = 
tv/u, the (complex) variable r travels (counterclockwise) full circle over the circular con- 
tour C, the diameter of which, of length 1/uj = T/(2ir), lies on the upper-half of the 
complex r-plane, with its lower end at the origin, r = 0, and its upper end at r = i/uo. 
Hence if the solution £(r) of (j2.5|) which emerges from some assigned initial data C(0), C'(0) 
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is holomorphic, as a (./V-vector-valued) function of the complex variable r, in the closed 
circular disk C encircled by the circle C in the complex r-plane, then the corresponding 
solution z(t) of (|2.1a|) . related to £(r) by (2.4), is completely periodic in t with period T, 
z(t + T) = z(t) (see (|2.2|) : and note that C( r )> considered as function of the real varia- 
ble t, is then periodic with period T/2, but z{t) is instead antiperiodic with period T/2, 
z(t + T/2) = —z(t), due to the prefactor exp(— iut), see (|2.4aj0 . 

In [1] it was proven that there indeed exists a domain, having infinite measure in 
phase space, of initial data C(0)j C'(0) such that the corresponding solutions ((t) of ()2.5|) 
are holomorphic in r in the disk C — and this fact implies the existence of an open 
domain, having as well infinite measure in phase space, of corresponding initial data z(0), 
i(0) such that the corresponding solutions z(t) of (|2.1aj) are completely periodic with 
period T, z{t + T) = z(t). In the following Section 3 we show that the singularities of 
the solutions ((t) of H2.5|) — considered as functions of the complex variable r — are 
branch points of square-root type, and in Section 4 we infer from this that, whenever the 
solution £(r) of (|2.5|) has a finite number of such branch points inside the circle C - 
generally nested inside each other, namely occurring on different Riemann sheets — then 
the corresponding solution z{t) of ()2.1aj) . considered as a function of the real "time" 
variable t, is again completely periodic, albeit now with a period which is an integer 
multiple of T. We also infer that when instead the solution £(r) of ()2.5|) has an infinite 
number of such square-root branch points inside the circle C — again, generally nested 
inside each other, namely occurring on different Riemann sheets — then the corresponding 
solution z{t) of (|2.1a|) . considered again as a function of the real "time" variable t, may 
be not periodic at all indeed it generally behaves chaotically (actually, strictly speaking, 
the chaotic character is not a property of a single solution, it rather has to do with the 
difference among the long-time behavior of a solution and those of other solutions which 
emerge from almost identical initial data — see below). 

Finally let us report a remarkable property of the system of ODEs (|2.1a|) . First of all 
we recall a trivial result, namely that the center of mass of the system (|2.1aj) . 

TV 

z(t) = N- 1 Y, z n(t), (2.7a) 

n=l 

rotates uniformly with period T, see (|2.2(l . since clearly these equations of motion, (|2.1a|) . 
entail 

iE + u 2 z = 0, (2.7b) 

hence 

z(t) = z(0) cos(wi) + l(0)(o;) _1 sin(wi). (2.7c) 
A less trivial result is that the sum of the squares of the particle coordinates, 

N 

zW(t)=N- 1 J2izn(t)} 2 , (2.8a) 

n=l 

evolves as well periodically, with period T/2. Indeed, as shown below, the equations of 
motion (J2.1a|) entail 

z& + (2uj) 2 zW = (4/N)H, (2.8b) 
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where the Hamiltonian H, see (|2.1b|) . is of course a constant of motion. Hence 

z^(t) = aexp(2zwi) + b + cexp(-2iut), (2.8c) 
where the three constants a, b, c are of course related to the initial data as follows: 



N 



r(0) = a + b + c = N~ l ^[^(0)] 2 



(2.9a) 



n=l 



N 



1(0) = 2iuj{a -b) = (2/JV) ^ i n (0)z n (0), (2.9b) 

71=1 

f( 2 )(0) = -(2w) 2 (a + b) = -(2lu) 2 z (2 \o) + (A/N)H. (2.9c) 

There remains to prove that the equations of motion Q2.1a[) entail (|2.8b|) . Indeed by 



differentiating twice the definition (|2.8aj) one gets 



N 



2 

N 



(2.10a) 



71=1 

and by using the equations of motion Q2.1a|) this yields 

JV N 



|(2) 



2 

N 



^ ^ (^"n ^ ^7i ) "I - ^ ^ 9nmZn(.Zn %n 
n,m=l; m^n 



\-3 



n=l 



(2.10b) 



Hence (using again the definition (|2.8a|) . and the symmetry of the coupling constants, see 



^) + ^ 2 z^ = jj 



N N 



\-2 



n=l 



n,m=l; m^n 



(2.10c) 



and using again (|1.2b[) . as well as the Hamiltonian equations z n — p n entailed by (|2~Tb)l . 
it is immediately seen, via (|2.1b|) . that the right-hand side of (|2.10cj) coincides with the 
right-hand side of ()2.8b|) . But the left-hand sides obviously coincide as well, so the result 
is proven. (We provided here a proof of this result for completeness, although of course 
this result can as well be directly inferred via the trick from the analogous finding recently 
proved for the system (2.5) [7]). 

3 The analytic structure of the solutions of ( 12.51 ) 

In this section we discuss the analytic structure of the solutions £(r) of the system of 
ODEs (|2.5j) . considered as functions of the complex variable r, and in particular we show 
that the singularities of these functions are branch points of square-root type (see below). 

These singularities are of course associated with values of the dependent variables Cn 
that cause the right-hand side of ()2.5|) to diverge, namely they are associated with "colli- 
sions" of two or more of the coordinates C, n . (By definition the "collision" of two or more 
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coordinates occurs when their values coincide; but we use inverted commas to underline 
that, since we are considering complex values of the independent variable r, such "colli- 
sions" need not correspond to actual collisions of the "particles" the motion of which as 
the real "time" variable t unfolds is described by the "physical" equations of motions <fl.l|> 
or (|2.1a|) — we shall return to this point below). Clearly the generic case corresponds to 
two-body collisions, occurring at some value r = r b such that, say, 

Ci(n) =Un). (3.ia) 

Here and below, without loss of generality, when discussing two-body collisions, we focus 
on the two particles carrying the labels 1 and 2. Note that it is natural to expect that 
this equation, 1)3. la)) , have one or more (possibly an infinity) of solutions T b in the complex 
r-plane for any generic solution £(r) of the ODEs 1)2.5)1 : while the equations characterizing 
a multiple collision, say 

Ci{n) = C2(n) = Csin) = ■■■ = Cuin), (3.1b) 

with 2 < M < N have generally no solution at all (for a generic solution £(t) of the 
ODEs (|2.5[) : there exist of course special solutions of the ODEs 1)2.5)1 which feature such 
multiple collisions, and we discuss them below to demonstrate that they as well feature 
branch points of square-root type; we moreover believe, due to the scaling character of the 
equations of motion 1)2.5)1 . that completely multiple collisions characterized by M = N are, 
as well as two-body collisions, featured by any generic solution of (|2.5|) — as suggested by 
the analysis of the three-body case, see Appendix B). 

To demonstrate that the singularity of the solution £(r) of the ODEs ()2.5|) associated 
with the collision ()3.1a)l is a branch point of square-root type we introduce the following 
ansatz, valid in the neighborhood of r = r b : 



Ci(r) = b + P(t - r h ) 1 ' 2 + v{r-T b ) + o(t - nf 2 + £ (r - n) 1 ' 2 , (3.2a) 

1=4 

00 

C 2 (r) =b-(3(r- n) 1 ' 2 + v(t - r b ) - a(r - r b f' 2 + £ of (r - r b ) 1 / 2 , (3.2b) 

1=4 

00 

Cn(r) =b n + v n (T-T b ) + ^2a\ n \r-T b ) 1 / 2 , n = 3, 4, . . . , N. (3.2c) 

1=4 

The fact that this ansatz is consistent with the two-body-collision condition 1)3. la|) is plain. 
It can moreover be shown (see Appendix A) that this ansatz, (3.2), is compatible with 
the ODEs (|2.5|) . and moreover that, while the coefficients j3 and Oij"^ featured by it are 
uniquely determined by the requirement that (3.2) satisfy 1)2.5)1 . the remaining coefficients, 
namely the 4 constants b, v, a and T b , see 

CE25), Cffl), and the 2(N - 2) constants b n , 
v n , see ()3.2c)) . are arbitrary (except for the obvious requirements b n ^ b, b n 7^ b m ) — and 
since altogether the number of these arbitrary constants is 2A^, we infer that this ansatz, 
(3.2), is adequate to represent in the neighborhood of r = r b the general solution of the 
system of second-order ODEs 1)2.5)1 . 

We therefore conclude that the singularities associated with the generic solution of 
()2.1a)) are branch points of square-root type, since clearly this is the singularity exhibited 
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by the right-hand sides of (3.2) — on the assumption that the infinite series they feature 
do indeed converge in a sufficiently small neighborhood of r = r&. 

It is easily seen that this same conclusion obtains if consideration is extended to include 
multiple collisions, see ()3.1b|) . Indeed the corresponding ansatz reads, in analogy to (3.2), 

Cn(r) = b + (3 n (r - n) 1 ' 2 + v(t - T b ) 

oo 

+ a n {T-T b f' 2 + ^a}" ) (r-T b ) / / 2 , n = l,2,...,M, (3.3a) 

1=4 

oo 

Cn(r) = b n + v n (T-T b ) + Y,a < i n) (T-n) 1 / 2 , n = M + l,M + 2,...,N. (3.3b) 

1=4 

The consistency of this ansatz with the ODEs 1)2.5)) is also demonstrated in Appendix A; 
but note that, as shown there, only the 3 constants b, v, r b appearing in the right-hand 
side of ()3.3a)l . and a common rescaling factor, say a, of the M coefficients a n , as well as 
the constants b n , v n appearing in the right-hand side of (|3.3b)l . are now arbitrary; so, the 
ansatz (3.3) features altogether 4 + 2(iV — M) = 2N — 2{M — 2) arbitrary constants. Hence 
for M > 2 this ansatz does not feature the full complement of 2N arbitrary constants 
required in order that (3.3) represent, for r ~ T b , the general solution of ()2.5)) — as indeed 
we expected (since we do not expect the generic solution of 1)2.5)1 to feature multiple 
collisions, at least with 2 < M < N). 

Let us end this section by reporting the similarity solution of 1)2.5)1 that indeed features 
an N-body collision. It reads, as can be easily verified, 

Cn(r) = b + p n (r- T h ) 1 ' 2 + v(T-n), (3.4a) 

with the 3 constants b, v and tj, arbitrary, while the iV constants (3 n are instead determined 
by the N algebraic equations 

N 

/3„ = -4 9n m (Pn-Pmy 3 . (3.4b) 

m=l, mj^n 



As it is well known (see for instance [2]), in the integrable equal-coupling-constants 



case, 



9nm = 9, (3.5a) 
the solution of ()3.4bl) is 

n = (-2 5 ) 1 / 4 e n , (3.5b) 
where the ./V real numbers £ n are the N zeros of the Hermite polynomial i5/jv(0 of degree N, 



H N (U) = 0. (3.5c) 

Note that the formula ()3.5b)) defines in fact 4 different sets of coefficients f3 n , due to 
the 4 possible determinations of the fourth root appearing in the right-hand side of this 
equation. 
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Finally let us emphasize that these findings entail that, even in the integrable equal- 
coupling-constants case, see (|3.5aj) . the solutions of the ODEs ()2.5|) do not possess the 
so-called Painleve property, namely they are not free of movable branch points, indeed 
their branch points (of square-root type, see for instance the special solution Q3.4a|) ) occur 
at values r = of the independent variable r which are not a priori predictable, but 
rather depend, in the context of the initial-value problem, on the initial data. But in 
the integrable case the number of these branch points is always finite: they are indeed 
determined by the requirement that the monic polynomial of degree N in the variable, 
say, C) the iV zeros of which are the iV coordinates Cn( r ) (and the coefficients of which are 
themselves polynomials of degree N in r), have a double zero (see for instance [2]). 

4 Periodicity of the solutions of the system of ODEs ( I2.1al ): 
theoretical considerations 

In this section we analyze the implications of the findings of the previous Section 3 as 
regards the periodicity of the solutions of the ODEs (|2.1a() . 

We know of course that, if a solution £(t) of (|2.5j) is holomorphic as a function of the 
complex variable r inside the circle C (see Section 2) — and we know that such solutions 
do exist, in fact in the context of the initial- value problem they emerge out of a set of 
initial data which has infinite measure in the phase space of such data [1] - then the 
solution z{t) of (|2.1aj) that corresponds to it via (2.4) is completely periodic with period T, 
see (I2~2l . 

z n (t + T) = z n (t). (4.1a) 

But the transformation (2.4) actually implies, as already noted in Section 2, an additional 
information, namely that in this case z{t) is completely antiperiodic with period T/2, 

z n (t + T/2) = -z n (t). (4.1b) 

Note that (|4.1b|) implies (|4.1a|) . while of course (|4.1aj) does not imply (|4.1b|) . 

Let us instead assume that a branch point of a solution ((t) of (|2.5|) . occurring, say, 
at r = Tb, does fall inside the circle C in the complex r-plane (see Section 2). Then, due 
to the square-root nature of this branch point, see (3.2), the evolution of the solution £(r) 
of 1)2. 5 Jl as the real time variable t unfolds is obtained by following the complex time- 
like variable r as it travels (2.4) along the circular contour C on a two-sheeted Riemann 
surface. Clearly the change of variable (|2.4a|) . (|2,4b|) entails then that the corresponding 
solution z{t) of ()2.1a|) is just as well completely periodic with period T, see (|4.1a|) . although 
in this case (|4.1bj) does no more hold. And of course this conclusion holds provided only 
one branch point of the solution £(t) of (|2.5|) falls inside the circle C in the main sheet 
of the Riemann surface associated with this solution, and no other branch point occurs 
inside the circle C in the second sheet of this Riemann surface, namely on the sheet entered 
through the cut associated with the branch point occurring inside C on the main sheet 
of the Riemann surface (of course this Riemann surface might feature many other sheets 
associated with other branch points occurring elsewhere hence not relevant to our present 
discussion) . 
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Let us now continue this analysis by considering, more generally, a solution £(r) of ((2.5(1 
that possibly contains more than one branch point inside the circle C in the main sheet of 
its Riemann surface (that do not cancel each other) so that by traveling along the circle C 
several additional Riemann sheets are accessed from the main sheet, and let us moreover 
assume that, on these additional sheets, additional branch points possibly occur inside 
the circle C which give access to other sheets, and that possibly on these other sheets 
there be additional branch points and so on. Let in conclusion B be the total number of 
additional sheets accessed by a point traveling around and around on the circle C on the 
Riemann surface associated with the solution £(t) of ((2.5(1 . This number B might coincide 
with the total number of branch points occurring, inside the circle C, on this Riemann 
surface — on all its sheets — or it might be smaller. Indeed, since each of these branch 
points is of square-root type, each of the associated cuts — if entered into — gives access 
to one additional sheet. But not all these sheets need be accessed; the total number B 
that are actually accessed depends on the structure of the Riemann surface, for instance 
no additional sheet at all is accessed if there is no branch point on the main sheet of the 
Riemann surface — even though other branch points may be present inside the circle C on 
other sheets of the Riemann surface associated with the solution Q(t) of ((2.5(1 . (It might 
also be possible that different branch points cancel each other pairwise as is the case for 
two branch points that are on the same sheet inside the circle C and generate a cut that 
starts at one of them and ends at the other). In any case the overall time requested 
for the point r(t) traveling on the Riemann surface to return to its point of departure 
(say, t(0) = on the main sheet) is [B + l)T/2, since a half-period T/2, see (|2.4b|) . is 
required to complete a tour around the circle C on each sheet, and the number of sheets 
to be traveled before getting back to the point of departure is overall B + 1 (including 
the main sheet). Hence in this case, as the real time t evolves, the solution £(t) of ((2.5(1 
will be completely periodic with period (B + l)T/2. Hence (see (2.4)) if B is even the 
corresponding solution z_(t) of ((2.1a() will be completely antiperiodic with the same period 
(B + l)T/2, z[t + (B + l)T/2] = —z(t), and completely periodic with the "odd" period 
(B + 1)T, z[t + (B + 1)T] = z(t). If instead B is odd, the solution ((t) of ((2~5|) as well 
as the corresponding solution z(t) of 1(2. la|) will both be completely periodic in t with the 
period (B + l)T/2 (which might be "even" or "odd" — of course, as an integer multiple 
of the basic period T), z[t + (B + l)T/2] = z(t) (so, in this case, the trajectories of z(t) 
will display no symmetry, in contrast to the previous case). 

In this analysis the assumption was implicitly understood that the total number B of 
additional sheets accessed by traveling around and around on the circle C on the Riemann 
surface associated with the solution £(r) of ((2.5(1 be finite (a number B which, as we just 
explained, might coincide with, or be smaller than, the total number of branch points of 
that Riemann surface that are located inside the circle C in the complex r-plane); and 
moreover we implicitly assumed that no branch point occur exactly on the circle C. Let 
us now elaborate on these two points. 

If a branch point occurs exactly on the circle C, then the "physical" equations of 
motion ((2.1a() become singular, due to a particle collision occurring at the real time t c 
defined mod (T/2) (see ((2.2(1 ) by the formula 

r b = [exp(2ioot c ) - l]/(2iw). (4.2) 

Indeed it is easy to check via (2.4) that the condition that t c be real coincides with the 
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requirement that the corresponding value of t&, as given by (|4.2jl . fall just on the circular 
contour C in the complex r-plane. The singularity is of course due to the divergence, at 
the collision time t = t c , of the right-hand side of the equations of motion (|2,la|k there is 
however no corresponding divergence of the solution z(t), which rather has a branch point 
of square root type at t = t c , see (3.2). But of course this entails that the speeds of the 
colliding particles diverge at the collision time t = t c proportionally to \t — i c | _1//2 , and 
their accelerations diverge proportionally to \t — t c |~ 3 / 2 . 

There is no a priori guarantee that the number of branch points inside C of a solu- 
tion £(t) of 1)2. 5 j) be finite, nor that the number B of additional sheets accessed according 
to the mechanism described above by moving around the circle C on the Riemann surface 
associated with that solution £(r) of ()2.5|) be finite (of course B might be infinite only if 
the number of branch points inside C is itself infinite). Obviously in such a case (B = oo), 
although the complex number r = r(t), see (|2.4b|) . considered as a function of the real 
"time" variable t, is still periodic with period T/2, neither the solution £(t) of (|2.5|h nor 
the corresponding solution z{t) of (|2.1aj) . will be periodic. The question that might then 
be raised is whether such a solution — in particular, such a solution z(t) of the "physical" 
Newtonian equations of motion (jl.lj) corresponding to the many-body problem character- 
ized by the Hamiltonian (1.2) — displays a "chaotic" behavior, namely, in the context of 
the initial-value problem, a "sensitive dependence" on the initial data. We shall return to 
this question below. 

So far we have discussed the relation among the analytic structure of a solution £(r) 
of (|2.5|) and the corresponding solution z(t) of (|2.1aj) . Let us now return to the simpler cases 
considered at the very beginning of this analysis and let us consider how the transition 
from one of the two regimes described there to the other occurs in the context of the 
initial- value problem for (J2.1a|) . and correspondingly for ()2.5|) . see (2.6). Hence let us 
assume again that the initial data for (|2.1a|) . and correspondingly for (|2.5|) (see (2.6)), 
entail that no branch point of the corresponding solution £(r) of (|2.5|) occurs inside the 
circular contour C on the main sheet of the associated Riemann surface, so that the 
corresponding solution z(t) of ()2.1a|) satisfies both (|4.1a() and (|4.1b|) . Let us imagine 
then to modify with continuity the initial data, for instance by letting them depend on an 
appropriate scaling parameter (a particular way to do so will be introduced in the following 
Section 5, as a convenient technique to present numerical results). As a consequence the 
branch points of the solution £(t) of (|2.5|) move, and the Riemann surface associated to 
this solution £(r) of 1)2. 5 Jl gets accordingly modified. We are interested in a movement of 
the branch points which takes the closest one of them on the main sheet of the Riemann 
surface from outside to inside the circle C. In the process that branch point will cross the 
circle C, and the particular set of initial data that correspond to this happening is then 
just a set of initial data that entails the occurrence of a collision in the time evolution of the 
many-body problem (|2.1a|) . occurring at a real time t = t c defined by (|4.2|) . as discussed 
above. After the branch point has crossed the contour C and has thereby entered inside 
the circular disk C, the corresponding solution z(t) of (|2.1a|) is again collision-free but its 
periodicity properties are changed. One might expect that the new solution continue then 
to satisfy 1)4. la)) but cease to satisfy (|4.1b|) . This is indeed a possibility, but it is not the 
only one. Indeed, since the time evolution of the solution z(t) of 1)2. la|) obtains via (2.4) 
by following the time evolution of the corresponding solution £(t) of (|2.5|) as the point 
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t = r(t) goes round the circle C on the Riemann surface associated with that solution, the 
occurrence of a branch point inside the circle C on the main sheet of that Riemann surface 
entails that the access is now open to a second sheet, and then possibly to other sheets if, on 
that second sheet, there also are branch points inside the circle C. If this latter possibility 
does not occur, namely if on that second sheet there are no branch points inside the 
circle C, then indeed there occurs for the corresponding solution z(t) of 1)2. la)) a transition 
from a periodicity property characterized by the validity of both 1)4. la|) and (|4.1h)l . to one 
characterized by the validity of 1)4.1 a|) but not of l|4.1hj) . If instead there is a least one 
branch point in the second sheet inside the circle C, then the periodicity — if any - 
featured after the transition by the solution z(t) of ()2.1aj) depends, as discussed above, on 
the number B of sheets that are sequentially accessed before returning — if ever — to the 
main sheet. 

To simplify our presentation we have discussion above the transition process by taking 
as point of departure for the analysis the basic periodic solution — that characterized 
by the validity of both ()4.1a)) and l)4.1b)l . the existence of which has been demonstrated 
in [1] — and by discussing how a continuous modification of the initial data may cause 
a transition to a different regime of periodicity, with the transition occurring in correspon- 
dence to the special set of initial data that yields a solution characterized by a particle 
collision, namely a set of initial data for which the Newtonian equations of motion become 
singular at a finite real time t c (defined mod(T/2)). But it is clear that exactly the same 
mechanism accounts for every transition that occurs from a solution z{f) of 1)2. la|) char- 
acterized by a type of periodicity to a solution z(t) of ()2.1a)) characterized by a different 
periodicity regime — or by a lack of periodicity. 

The final point to be discussed is the question we postponed above, namely the char- 
acter of the nonperiodic solutions z(t) of ()2.1a)l (if any), which we now understand to 
be characterized, in the context of the mechanism described above, by access to an end- 
less sequence of different sheets — all of them generated by branch points of square-root 
type — of the Riemann surface associated with the corresponding solution £(r) of 1)2.5)1 . 
The following two possibilities can be imagined in this connection — which of course does 
not entail they are indeed both realized. 

The first possibility — which we surmise to be the most likely one to be actually 
realized — is that an infinity of such relevant branch points occur quite closely to the 
circular contour C, hence that there be some of them that occur arbitrarily close to C. 
This then entails that the corresponding nonperiodic solutions z(t) of 1)2 .la)) manifest 
a sensitive dependence on their initial data (which we consider to be the signature of 
a chaotic behavior). Indeed a modification, however small, of such initial data entails a 
modification of the pattern of such branch points, which shall cause some of them to cross 
over from one to the other side of the circular contour C . But then the two solutions z(t) 
of ()2.1a)) corresponding to these two assignments of initial data — before and after the 
modification, however close these data are to each other — will eventually evolve quite 
differently, since their time evolutions are determined by access to two different sequences 
of sheets of the Riemann surfaces associated with the two corresponding solutions C(r) 
of 1)2 .5|) — two Riemann surfaces which themselves need not be very different (to the 
extent one can make such statements when comparing two objects having as complicated 
a structure as a Riemann surface with an infinite number of sheets produced by an infinite 
number of branch points of square-root type). So this is the mechanism whereby a chaotic 



496 



F Calogero and M Sommacal 



behavior may develop for the system Ql.ljl — but of course not in the integrable case with 
equal coupling constants (in that mentioned above, the number of branch points 

of the solutions £(r) of (j2,5j) is always finite, since the N coordinates Cn( r ) are i n this case 
the N zeros of a polynomial of degree N the coefficients of which are polynomials in the 
variable r [2]). Note that the emergence of such a chaotic behavior would not be associated 
with a local exponential divergence of trajectories in phase space — it would be rather 
analogous to the mechanism that causes a chaotic behavior in the case of, say, a triangular 
billiard with angles which are irrational fractions of tt — a chaotic behavior also not due 
to a local separation of trajectories in phase space, but rather to the eventual emergence 
of a different pattern of reflections (indeed of any two such billiard trajectories, however 
close their initial data, one shall eventually miss a reflection near a corner which the other 
one does take, and from that moment their time evolutions become quite different). 

A different possibility, which we consider unlikely but we cannot a priori exclude at this 
stage of our analysis, is that nonperiodic solutions z(t) of (|2.1aj) exist which are associated 
with a Riemann surface of the corresponding solutions £(t) of 1)2. 5 Jl that, even though 
it possesses an infinite number of relevant branch points inside the circular contour C, 
it features all of them — or at least most of them, except possibly for a finite number 
of them — located in a region well inside C, namely separated from it by an annulus 
of finite thickness. Clearly in such a case two nonperiodic solutions z(t) of 1)2. la)) which 
emerge from sufficiently close initial data separate slowly and gradually throughout their 
time evolution, hence they do not display a sensitive dependence on their initial data — 
hence, in such a case there would be solutions which are nonperiodic (nor, of course, 
multiply periodic) but which nevertheless do not display a chaotic behavior — or, to be 
more precise, there would be sets of initial data, having nonvanishing measure in the phase 
space of initial data, which yield such nonperiodic (yet nonchaotic) solutions. 

As indicated by its title, this paper is mainly focussed on the periodic solutions. In 
the following Section 5 the analysis of their phenomenology given in this section is com- 
plemented by the display of numerical simulations. We shall also exhibit solutions which 
appear nonperiodic and perhaps chaotic, although such characteristics can of course never 
be demonstrated with complete cogency via numerical examples. 



5 Periodicity of the solutions of the system of ODEs ( I2.1al ): 
numerical simulations 

In this section we display several numerical solutions of the Newtonian equations of motion 
(f2~Ta)l with 

u) = 2tt (5.1a) 
hence (see JH2J)) 

T = 1. (5.1b) 

These results confirm graphically the findings discussed in the preceding Section 4. The 
strategy of our presentation is to exhibit a sequence of numerically-computed solutions of 
(|2.1a|) with (|5.1aj) corresponding to different choices of initial data, for various models 
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characterized by an assignment of the number N of particles and of the values of the 
coupling constants g nm (that are always assumed to satisfy the symmetry property (|1.2b|) , 
9nm = 9mn)- For obvious reasons of simplicity, see below, we restrict consideration to 
N = 3 and N = 4. For each model we consider sequences of motions characterized by sets 
of initial data linked to each other by the formulas 

x n (0) = A- 1 ^ ), y n (0) = A- 1 ^, (5.2a) 
x n (0) = A4°)-2^(A-A- 1 )yi°), y n (0) = \y® + 2tt(A - A^xf , (5-2b) 

of course with z n = x n + iy n . Here A is a positive rescaling parameter the different values 
of which identify different sets of initial data (while the data Xn , Vn^ , Xn^ , Vn^ are kept 
fixed). The motivation for such a choice is that, as it can be easily verified by rescaling 
appropriately the dependent and independent variables (see (2.3) and (2.6) with (|5.1a|l ). 
these sets of initial data (5.2) identify different solutions z(t) of the system of ODEs (|2.1aj) 
with (|5.1a|) that correspond to different solutions £(r) of the system of ODEs (|2,5|) related 
to each other by the following change of initial data, 

£(0) = A" 1 ^, C'(0) = AC' (0) . (5.2c) 

This, as can be easily verified, entails these solutions £(r) are related to each other merely 
via a rescaling of dependent and independent variables by a constant factor. Hence all 
these solutions are associated to the same Riemann surface except for a shrinking of the 
complex r-plane by a common factor — which, as can be readily verified, is just A 2 . 
Therefore, in the context of the discussion of the preceding Section 4, to analyze the 
motions yielded by the initial data (5.2) one can just imagine to multiply the diameter 
of the circle C by the factor A 2 without modifying the Riemann surface — so that larger 
values of A entail that more branch points fall within C. 

So for small enough values of A — namely for initial conditions characterized by large 
particle coordinates hence by large interparticle separation and by large initial velocities 
almost "orthogonal" to the initial positions (see (5.2a,b)) — the circle C shrinks to a 
small enough radius so that it contains no branch points, hence the corresponding motion 
(solution of (|2.1aj) with (5.2a,b)) is completely periodic with period T = 1, see (|4.1aj) 
and (5.1b), and moreover it features the symmetry property (|4.1b|) . As A gets increased 
one goes through the scenarios described in Section 4 — since, as we concluded above, 
such an increase can be interpreted as amounting merely to an increase of the radius of 
the circle C, causing thereby more and more branch points to be enclosed by it. For very 
large A this process can cause all branch points to be enclosed inside C, a situation which 
is of course equivalent, as regards the periodicity of the solutions of as functions of 
the real time variable t, to all of them being outside — hence the expectation, that for very 
large A the motion be again completely periodic with period T = 1 and also symmetrical, 
see (|4.1ajl . (|4.1b|) (but this need not necessarily happen, since the possibility that the 
solution £(t) possesses branch points at infinity cannot be excluded). This expectation is 
indeed confirmed by some of the examples reported below, while in other cases for very 
large A the numerical simulations become too difficult to be performed reliably, and the 
corresponding trajectories become unsuitable for transparent display — note that for very 
large A the particles are initially all very close to the origin (hence to each other) and with 
large initial velocities, see (5.2a,b). 
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But in this paper we do not elaborate on the numerical aspects, except to reassure the 
skeptic reader that we made sure in each case of the reliability of the results presented 
below; we refer for more details to [6]. Let us however note that, while one might think 
that a useful check of the accuracy of the computation is provided by a verification that the 
numerically evaluated coordinates z n (t) indeed imply that the collective coordinates z(t) 
respectively z( 2 )(i), see ()2.7a|) respectively (|'2.8a(l . evolve periodically with periods T = 1 
respectively T/2 = 1/2 according to (|2.7c|) respectively (|2.8cj) . these tests are in fact not 
at all cogent: indeed, even poorly evaluated coordinates z n (t) tend to yield collective 
coordinates z{t) and z^ (t) that evolve properly. This happens because the equations of 
motion that determine the evolution of the collective coordinates z(t) and z^ 2 \t) are in 
fact so simple (see (2.7) and (2.8)), that numerical errors made in the integration of (|2.1al) 
tend to cancel out when z(t) and z^ (t) are evaluated. 

Let us also note that, in all the examples considered below, we assigned for simplicity 
real integer values to the coupling constants; but we did check in every case that the qual- 
itative character of the motions does not change if these coupling constants are replaced 
by (neighboring) values which are neither entire nor real. 

The first example we consider is characterized by the following parameters: 



N = 3; gu = g 2 i = 1, 523 = 932 = 10, 531 = 913 = -2, (5.3a) 

)y the folk 
data via (5.2): 



and by the following values of the parameters Xn^ , Un^ , xffl , ijn^ characterizing the initial 



T (°) 

x l 


= 1, 


(0) 


0, 


±i 0) = 0, 


T (0) 

x l 


= 0, 


(0) 


1, 


*?> = -2 


T (0) 

x l 


= 0.5, 


(0) 


= 0.5, 


if = 



y? = -1, 
yf ] = 0, 

y? = 0. (5.3b) 

The following Table 5.1 provides an overview of the main features (periodicity, sym- 
metry) of the motions that emerge from initial data determined by (5.2) with (|5.2a|) . as 
well as an indication of the figures (if any) where the corresponding trajectories are dis- 
played. Here and always below the trajectories of particle 1, 2 respectively 3 are shown 
in red, green, respectively blue (and below, those of particle 4 in yellow). Whenever we 
felt such an additional indication might be usefully displayed we indicated with a black 
diamond the initial position of each particle (at t = 0), with a black dot the position at 
a subsequent time t = t\ (generally chosen to coincide with some fraction of the period, 
t\ = T/p, with p a conveniently chosen positive integer), and with smaller black dots the 
position at every subsequent integer multiple of t± (namely at t = i& = kt\, k = 2, 3, . . .); 
in this manner the direction of the motion along the trajectories can be inferred (from the 
relative positions along the trajectories of the diamond and the larger dot), as well as some 
indication of the positions of the particles over time, as they move (by counting the dots 
along the trajectory). Of course a much more satisfactory visualization of the motions is 
provided by simulations in which the particle motions are displayed as they unfold over 
time (as in a movie); it is planned to make available soon, via the web, the numerical code 
suitable to perform such simulations on personal computers [6]. Let us emphasize that 
such simulations are particularly stunning to watch in the case of high-period trajectories, 
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which are very complicated (see below), so that the fact that the particles return eventu- 
ally exactly on their tracks appears quite miraculous and is indeed a remarkable proof of 
the reliability of the numerical computation. 



Table 5.1 
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Now some comments on these results (see Table 5.1 and the set of Figs. 5.1). For 
A < 0.712 there clearly are no singularities inside C, hence the motion is periodic with the 
basic period T = 1, and moreover antiperiodic with period T = 1/2, so that all trajectories 
are symmetrical (see for instance Figs. 5.1a,b). A transition occurs at some value of A larger 
than 0.712 but smaller than 0.81, and it causes a branch point to enter inside C, opening 
the way to a second sheet (of the Riemann surface associated with the solution £(t) of (2.5) 
with (5.3)) and thereby causing the trajectories of the solution z(t) of (2.1) to loose their 
symmetry, while still preserving period T = 1. This transition is due to a collision among 
particles 1 (red) and 3 (blue), occurring at a time of the order of, or maybe a bit less than, 
(3/10)T = 3/10 (see Figs. 5.1b and 5.1c, where clearly h = T/10 = 1/10). Another, more 
dramatic, transition occurs for a value of A between 0.9 and 0.97, and it causes a major 
increase in the complication of the motion, perhaps a transition to chaos. This transition 
is of course due to the entrance inside C of a branch point that opens the way — as the 
complex variable r travels around and around on the circle C — to a very large number 
of sheets of the Riemann surface associated with the solution £(r) of (2.5) with (5.3), 
possibly to an infinity of them, resulting in motions that look chaotic: this we indicate in 
Table 5.1 with the acronym HSL, which stands for "Hie Sunt Leones" — see the column 
with A = 0.97. This transition is due again to a collision between particles 1 (red) and 3 
(blue), at a time of the order, or maybe a little less, than (8/10)T = 8/10 (see Fig. 5. Id, 
where clearly again t\ = T/10 = 1/10). Further increases of A cause instead a decrease of 
the complication of the trajectories, characterized first of all by a return to periodic (if still 
very complicated, due to the large periods) motions, and subsequently by a progressive 
decrease of the period (see the relevant graphs, from Fig. 5.1e to Fig. 5.1m); this is of course 
interpreted as due to the fact that the access to additional sheets of the Riemann surface 
associated with the solution Q(t) of (2.5) with (5.3) gets shut off because also the second 
branch point corresponding to the cut that opened the way to those sheets gets enclosed 
inside the circle C. Thus the trajectories shown in Fig. 5.1e (corresponding to A = 0.974; 
see Table 5.1) are interpreted, on the basis of the discussion of the preceding Section 4, 
as corresponding to access to 30 sheets altogether of the Riemann surface associated with 
the solution £(r) of (2.5) with (5.3); likewise those shown in Fig. 5. If (corresponding to 
A = 0.99; see Table 5.1) are interpreted as corresponding to access to 22 sheets of the 
Riemann surface associated with the solution £(r) of (2.5) with (5.3); those shown in 
Fig. 5.1k (corresponding to A = 1.2; see Table 5.1) are interpreted as corresponding to 
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access to 2 sheets, and finally those shown in Fig. 5.1m (corresponding to A = 5; see 
Table 5.1) are interpreted as corresponding to access to just 1 sheet, namely just the main 
one. 

The second example we consider is characterized by the following parameters: 

N = 3; 512 = 521 = 10, 523 = 532 = 3, 531 = 513 = -10, (5.4a) 

and by the following values of the parameters xffl , Vn^ , Xn^ , characterizing the initial 
data via (5.2): 
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Since the qualitative picture is analogous to that displayed above, we limit our presen- 
tation to the display of Table 5.2 and of the corresponding set of graphs (see Figs. 5.2), 
letting to the alert reader the fun to repeat the discussion as given above. We merely 
emphasize that, in contrast to the previous case, we see in this case also trajectories with 
an even period, as well as some high-period trajectories that are symmetrical (the alert 
reader will have no difficulty in identifying the number of sheets of the Riemann surface 
associated with the solution ("(r) of (2.5) with (5.4), for each of the cases reported in 
Table 5.2). One is also led to suspect, from the data reported in the last columns of 
Table 5.2, that the Riemann surface associated with the solution £(r) of (2.5) with (5.4) 
possesses (at least) 7 branch points at or near infinity (of course on different sheets, so 
that they do not cancel each other). 

The third example we consider is a 4-particle one, characterized by the following pa- 
rameters: 



N = 4; 5 i2 : 

523 = 532 = -9, 



921 = -3, 513 

524 = 542 = 12, 



531 = 9, 514 
534 = 543 = 3, 



541 



-12, 



(5.5a) 



and by the following values of the parameters Xn^ , yn^ , xffl , ijn^ characterizing the initial 
data via (5.2): 
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A representative selection of our findings is reported in Table 5.3, and some of the corre- 
sponding trajectories are displayed in the Figures as indicated there. 
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For A < 0.533 we have no singularities inside C, hence the motion is periodic with the 
basic period T = 1, and moreover antiperiodic with period T = 1/2, so that all trajectories 
are symmetrical (see Figs. 5.3a,b). A transition occurs at some value of A between 0.533 
and 0.8, due to a collision among particles 3 (blue) and 4 (yellow), occurring at a time 
of the order of, or maybe a bit less than, (3/10)T = 3/10 (see Figs. 5.3b and 5.3c, where 
clearly t\ = T/10 = 1/10). Another transition occurs for a value of A greater than 0.9687 
but less than 0.98, and it causes a branch point to enter inside C at some value of A 
between 0.8 and 0.9687, opening thereby the way to a third sheet of the Riemann surface 
associated with the solution £(t) of (2.5) with (5.5). This second transition is due to 
a collision, in this case among particles 2 (green) and 4 (yellow), occurring at a time 
of the order, or maybe a little less, than (8/10)T = 8/10 (see Fig. 5. Id, where clearly 
again t\ = T/10 = 1/10). Then two more transitions occur (see Fig. 5.3e,f,g) — one 
for A between 0.98 and 0.985, the other for A between 0.985 and 0.99 — each opening the 
way to one more sheet of the Riemann surface associated with the solution £(r) of (2.5) 
with (5.5) (recall the discussion of Section 4: overall access to P sheets including the main 
one, P = 1 + B in the notation of Section 4, yields a completely periodic symmetrical 
motion with period P if P is odd — with every sheet being visited twice in each period] 
it yields a completely periodic unsymmetrical motion with period P/2 if P is even - 
with every sheet being visited once in each period: hence for the trajectories shown in 
Figs. 5.3a,b P = 1, for those shown in Figs. 5.3c,d P = 2, for those shown in Fig. 5.3e 
P = 3, for those shown in Fig. 5.3f P = 4, for those shown in Fig. 5.3g P = 5). A more 
dramatic transition occurs for a value of A between 1 and 1.01, and it causes a major 
increase in the complication of the motion, perhaps a transition to chaos (see Fig. 5.3h). 
Like in the first example, and of course due to the same mechanism as discussed above, 
further increases of A (A > 1.225) cause instead a decrease of the complication of the 
trajectories, characterized first of all by a return to periodic motions, and subsequently 
by a progressive decrease of the period (see Figs. 5.3i-m). 

The fourth and last example we report is characterized by the following parameters: 

N = A- 512 = 521 = 5, 5i3 = 531 = 10, 514 = 541 = 20, 

523 = 532 = -20, 524 = 542 = -10, 534 = 543 = -5, (5.6a) 

and by the following values of the parameters xffl , y^ , x$ , characterizing the initial 
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data via (5.2): 
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A representative selection of our findings is reported in Table 5.4, and some of the corre- 
sponding trajectories are displayed in the Figs. 5.4 as indicated there. 



Table 5.4 



A 


0.5 


0.736 


0.8 


0.894 


0.9 


0.904 


0.905 


0.91 


0.95 


Period 


1 


1 


1 


1 


2 


2 


HSL 


HSL 


HSL 


Symmetry 


Yes 


Yes 


No 


No 


No 


No 








Fig. 5.4 


a 


b 


c 


d 


e 












A 


1 


1.1 


1.125 


1.13 


1.15 


1.2 


1.5 


2 


5 


10 


Period 


HSL 


HSL 


2 


2 


1 


1 


1 


1 


1 


1 


Symmetry 






No 


No 


No 


No 


No 


No 


No 


No 


Fig. 5.4 


f 






g 




h 






i 





Since the qualitative picture is analogous to those displayed above, we limit our pre- 
sentation to the display of Table 5.4 and of the corresponding set of Figs. 5.4, letting to 
the alert reader the fun to repeat the discussion as given above. 



6 Outlook 

The results proven in this paper provide an additional explicit instance of a phenomenon 
whose rather general scope has been already advertized via a number of other examples, 
treated elsewhere in more or less complete detail [1-4]. 

An analogous treatment to that given in this paper will be published soon [5] in the con- 
text of the "generalized goldfish model" [8, 2], which is somewhat richer inasmuch as it 
features branch points the nature of which depends on the values of the coupling constants 
(in contrast to the case treated herein, where all relevant branch points are of square-root 
type), and has moreover the advantage that its treatment in the CV context is directly 
interpretable as a genuine (i.e., rotation- invariant) real many-body problem in the plane 
(as discussed at the end of the introductory Section 1, this is not the case for the model 
considered herein — although there does exist a modified version of it in which this 
"defect" is eliminated [9]). On the other hand the class of many-body problems with 
inverse-cube interparticle potentials, as considered herein, have been (especially, of course, 
in the integrable version with equal coupling constants) much studied over the last quarter 
century, while the "goldfish" model has not yet quite acquired a comparable "classical" 
status. 

Of course a more detailed, if perhaps less "physical" , understanding of the dynamics of 
the model studied herein might be gained by investigating numerically the solutions ((t) 
of the equations of motion (2.5) as functions of the complex variable r (rather than the 



Periodic Solutions of a System of Complex ODEs. II. Higher Periods 



509 



solutions z{t) as functions of the real variable t), and in particular by mapping out the 
detailed shape of the multi-sheeted Riemann surfaces associated with these solutions; and 
analogous considerations apply to all the models [1-5] in which the "trick", see Section 2, 
plays a key role. This remains as a task for the future. 
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Appendix A: the ansatz (3.2) 

In this appendix we show that the ansatz (3.2) is consistent with the ODEs 1)2.5(1 . and that 
the requirement that (3.2) satisfy (|2.5[) allows to determine uniquely the coefficient (3 and 

(n) 

in principle as well all the coefficients a\ featured by this ansatz, while the remaining 
coefficients, namely the 4 constants b, v, a and in ()3.2a|) . ()3.2b|) . and the 2(N — 2) 
constants b n , v n in (|3.2cj) . remain arbitrary (except for the obvious requirements b n ^ b, 
b n 7^ b m , which are hereafter assumed to hold). 

Indeed the insertion of (3.2) in (|2.5j) with n = 1 respectively n = 2 yields 



1=4 




3 



+ Yl dlmib-b^-Hl + ib-bm)- 1 f3(T-T b ) 




v - v m )(r - n) 




(A.la) 



respectively 




nr 3/2 - \a{r - r fe )- 1/2 + E " 



1=4 




-3 



+ 




-1 




v - v m )(T - n) 
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a(r - nf' 2 + £ aj 



Z=4 



„( 2 ) 



a 



(m) 



(t-t 6 



.J/2 



(A.lb) 



We now equate in each of these two equations the coefficients of the term (r — t&) 3 / 2 
and we get the (same) relation that determines the coefficient /3, 

= (W2) 1/4 . (A.2) 

And we also see that, with this assignment of (3, in both equations, (|A.la|) . (|A.lb|) . the 
terms that multiply (r — T b )~ 1 / 2 also match exactly. 

It is also easily seen that the remaining terms can be matched recursively (by expanding 
the right-hand sides of these equations, l|A.la|) . (|A.lb|) — as well as those of (|A.lc|) . see 
below — in powers of (r — r^) 1 / 2 ), and that in this manner the coefficients ap get 
uniquely determined, for instance 



a f = Yl i( 13 9jm ~ 3<7 i+ i jm )/20](6 - b 



N 



3 

m ) i 



j = 1,2 mod(2). 



(A.3a) 



m=3 



To complete the analysis one must also check the remaining equations (|2.5|) . with n > 2. 
Insertion of the ansatz (3.2) in these yields 



^M a («) (r _ T6) (H/2 



1=4 



gm(b n -b)-ni + {b n -b)- 1 



P(t - n) 1 / 2 + (v n - v)(r - n ) 



1=4 



+ a(r-r b f/ 2 + Y,{a ( r ) -a < i 1) )(r-r b ) 



a/2 



+ gn2{b n -b)-H l + (b n -b) 



-i 



(3{t - Th ) 1 ' 2 + {v n - V){T - T h ) 

-3 



a/2 



1=4 



N 



+ E 9nm(b n -b m ) 3 <l + (b 

m=l, m^n V 



n u m J 



(V n ~ V m )(T - Tb) 



1/2 



1=4 



n = 3, 4, . . . , N. 



(A.lc) 



The consistency of these equations is plain, as well as the possibility they entail to compute 
in principle — in conjunction with (|A.lb|) . l)A.lc|) — the coefficients a\ ■ For instance 
one immediately obtains 



a 



(n) 



N 



(gnl+9n2)(b n -b) 3 + (bn ~ K 

m=3, m^n 

n = 3,4,... , N. 



\-3 



(A.3b) 
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The treatment of the multiple-collision case, see (|3.1b|) . via the corresponding ansatz, 
(3.3), is closely analogous, but in this case the matching of the coefficients of the terms 
proportional to (r — 17,) ~ 3 / 2 respectively to (r — it,) -1 / 2 yields 

M 

0n = -4 Yj 3nm(Pn ~ Pm)~ 3 , U = 1, 2, . . . , M, (A.4a) 

m=l, m^n 
M 

a n = -4 ^2 g nm (Pn- Pm)~ A {a n ~ a m ), n = l,2,...,M. (A.4b) 

m=X, mj^n 

So in this case the coefficients /3 n and a n get fixed by these algebraic equations (up, for 
the coefficients a n , to a common factor). But the basic conclusion about the square-root 
nature of the branch points is confirmed, as indeed entailed by the ansatz (3.3). 

And it can be easily verified that analogous conclusions also obtain in the case of more 
general multiple collisions in which the coordinates collide in groups, for instance 

(i(n) = bin), ( 3 (r b ) = ( 4 (n) = (s(n), ( n {n) ^ Ci(n), 

Cn(n) + Cafo), Cn(n) + Cm(n), n, m = 6, 7, . . . ,N, (A.5) 

the appropriate ansatz being in this case of type (|3.3a|) for the colliding coordinates Cn(T~), 
n = 1, 2, . . . , 5 (with the same values of the constants b and v for each subgroup of colliding 
coordinates, say, in self-evident notation, b\ = 62, &3 = &4 = &5, v\ = V2, v 3 = i>4 = fs), of 
type (|3~3bl for the non colliding coordinates Cn(i")> n = 6, 7, . . . , N. 

Appendix B: the three-body case 

In this appendix we discuss the solution of the evolution equations (|2.5|) in the three-body 
case, iV = 3. As we show below the solution of this three-body problem can be reduced 
to quadratures. This fact was discovered by C Jacobi almost two centuries ago (but this 
finding was only published after his death [10]), was then rediscovered by C Marchioro 
(who also did not publish this finding) and was finally neatly presented by D C Khandekar 
and S V Lawande [11], who actually treated (|2.1a|) rather than (|2.5|) (they were not aware 
at the time of the "trick" (2.4) that relates (|2.5|) to (|2.1a|) . and they mainly focussed on the 
integrable case with equal coupling constants and on its connection with the corresponding 
quantal problem — for a bit more on the history of this problem see Section 2.N of Ref. [2]). 
So we now focus on the evolution equations 

3 

Cn = ^ 9nm{Cn ~ Cm) ; fl = 1, 2, 3, Qnm = 9mm Cn = Cn( 7 ")) (B-1-) 

m=l, mj^n 

which clearly obtain from the Hamiltonian 

h = ~ ( P i + P i +p\) + ~ Ma - C2)- 2 + 923«2 - C3)- 2 +53i(c 3 - ar 2 ] 

hence entail that the quantity 

H = \ (Cf + + C 3 2 ) + \ [5l2(Cl " C 2 r 2 + 523(C2 - C 3 )" 2 + 531 (C3 - Cl)~ 2 ] (B.2b) 
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is a constant of the motion. It is also plain that the initial position, C(0), an d the (initial) 
velocity, V, of the center of mass 

C = (Ci + C2 + GOA (B.3) 



are two additional constants of motion, since the evolution equations (jB.lj) clearly entail 
that the center of mass moves uniformly: 

C(t) = C(0) + Vt. (B.4) 
And we moreover now know (see the last part of Section 2) that the quantity 

C (2) = (d 2 + Cl + CI) /3 (B.5a) 

evolves according to the simple equation (see (|2.8b|l ) 

(C (2) )" = (4/3)ff, (B.5b) 

with H defined by (IB~2bl . 

It is now convenient [11] to introduce the "Jacobi coordinates" 

r, = 2- 1 / 2 (Ci - C2), e = e-^CCi + C2 - 2Cs), (B.6a) 

and moreover to introduce the corresponding "circular coordinates" by setting 

rj = psm(6), £ = pcos(8). (B.6b) 
Note that these definitions entail 

Ci = C + 2- 1/2 r/ + = C - (2/3) 1 / 2 poos[0 + (2tt/3)], (B.7a) 

C2 = C " 2~ 1/2 r? + e- 1 /^ = C - (2/3) 1 / 2 pcos[0 - (2tt/3)], (B.7b) 

C 3 = C " (2/3) 1/2 £ = C " (2/3) 1 / 2 /0 cos 8, (B.7c) 
as well as 

H = (3/2)V 2 + (l/2)p' 2 + (l/2)pV 2 + (l/A)p- 2 f(9), (B.8a) 

f(0) = gu/ sm 2 (6) + g 23 / sin 2 [6 + (2vr/3)] + 531/ sin 2 [6 - (2tt/3)] , (B.8b) 

and 

C( 2 ) =C 2 + (l/3)/9 2 . (B.9) 
From the last formula and (|B.4|) . (|B,5b|) we get 

(p 2 )" = AH -6V 2 , (B.lOa) 

hence 

p 2 ( T ) = A (r 2 + 25r + C) = A(t - r+)(r - r_), (B.lOb) 

^ = 2F-3F 2 , (B.lOc) 

t± = -B±D, (B.lOd) 

D 2 = B 2 - C. (B.lOe) 
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Here we have introduced two other constants of integration, B and C, while the constants A 
and D are defined by (|B.10c|) and (|B.10ej) . 

We now introduce this expression of p(r) in (jB.8a|) and easily obtain the formula 

9> 2 = -[A 2 D 2 + (l/2)f(9)]p-\ (B.lla) 

that entails the quadrature 

j° dO' [—A 2 D 2 - (1/2) / ■(*')] ~ V2 = f dr' [p(r')] ~ 2 . (B.llb) 

The integral in the right-hand side of this formula, (|B,llb|) . is easily done (see (B.10)), 
and one arrives thus at the final formula 

F{9) = arccotan [(r + B)/D], (B.12a) 
F{9) = J d9'[l + (l/2)(AD)- 2 f(e')]- 1/2 , (B.12b) 

where of course f(9) is defined by (|B.8b|) , 

In the integrable "equal-coupling-constants" case, 

g nm = 9, (B.13a) 

it is easily seen that 

f(9) = 9 5 [sin(30)]- 2 , (B.13b) 
hence the integration in the right-hand side of (|B.12b|) is easily performed to yield 

F{9) = (l/3){£ - arcsin[cos(36')//i]}, (B.13c) 
where E is an integration constant and 

p=[l + (9/2)g(AD)- 2 ] 1/2 . (B.13d) 
Hence in this case one gets for 9{t) the completely explicit expression 

9 = (1/3) arccos [psm{E - 3 arccotan [(r + £)/£>]}], (B.14) 

which, together with ()B.10b|) and (B.4), provides via (B.7) completely explicit expressions 
of the coordinates Cn( T )j n = 1,2, 3. From these it is immediately seen that the coordina- 
tes CnC 7 ")) considered as functions of the complex variable r, feature two kinds of square-root 
branch points: those occurring at r = t±, see (|B.10b|) . (|B.10d|) . which clearly correspond 
to triple collisions, see (|B.10b|) and (B.7); and those occurring at r = t^, 

r b = -B + Z)cotan{(l/3)[ J B - arcsin(l/^) + (2vr/3)A:]}, k = 0,1,2, (B.15) 

which correspond instead to pair collisions. Indeed it is easily seen that, for r ~ 75, 

9 = (2ir/3)k + a[(T-n)/D] 1/2 + 0(\T-T b \), k = 0,1,2, (B.16a) 

a = i(2/3) 1 / 2 [1 + (2/9)(AD) 2 /g] ~ 1/4 D [(r b + B) 2 + D 2 ] ' 1/2 ; (B.16b) 
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and (|B.16a|) clearly entails, see (B.7), that a pair collisions occurs at r = tj. 

The analysis of the periodicity of the solutions of (2.1) in this case is plain; of course 
the fact that the number of branch points is finite entails that in this integrable case all 
nonsingular solutions are completely periodic, either with period T, see (|2.2(l . or with a 
period which is a (small) entire multiple of T, confirming the results implied by the general 
treatment [2] of the integrable equal-coupling-constant iV-body problem (2.1). 

In the general case with arbitrary coupling constants one can rewrite (IB.12hjl as follows 
(by setting u = 4 sin 2 6' — 3) : 

F{9) = (1/2) / duu(l -u)- 1 ' 2 



• u 2 {u + 3 + G) + G + G{3 + u) 3 ' 2 {l-u) l/2 \ 1/2 , (B.17a) 
G = {2g 12 - g 23 - g 31 ){AD)~ 2 , 

G = 9(g 2 3 + g 3 i)(AD)- 2 , G = 3 1 / 2 (g 2 3 - g 31 )(AD)- 2 , (B.17b) 

17(0) = 4[sin0] 2 -3. (B.17c) 

In the equal-coupling-constants case, see (|B.13a|) . things simplify because 

G = G = 0, (B.18a) 

see (|B.17b|) . hence the integral in the right-hand side of (|B.17|I can then be performed in 
terms of elementary functions, since via the change of variable 

u 2 (u + 3) = iv, (B.18b) 

somewhat miraculously (albeit not surprisingly) (B.17) yield 

F(0) = (l/3)/ dv {{1 - v){G + v)}- 1 ' 2 , (B.18c) 
W{6) = [sin(30)] 2 ; (B.18d) 

and of course one gets thereby again (|B.13c|) . 

If instead (only) the two coupling constants g 23 and g 3 \ coincide, 523 = 531; entailing 
G = 0, the integral in the right-hand side of (|B.17|) is of elliptic type (of course the same 
outcome, namely the elliptic character of the solutions [9], characterizes any case with 
two equal coupling constants, although this is not immediately evident from (jB.17|) ). But 
not much enlightenment seems to obtain, even in this (somewhat simpler) case, from an 
attempt to perform explicitly the integration in the right-hand side of ()B.17|) . 

One must rather try and evince information directly from the integral representation 
of F(9), see (|B.12b|) with (jB,8b|) . or (B.17), or from the following representation which 
the diligent reader will have no difficulty in deriving: 

/exp(2i0) 
dww- 1 (w 3 - 1) [(P(w)]- 1/2 , (B.19a) 

<p{w) = (w 3 - l) 2 - 4(AD)- 2 w{g 12 [w - exp(2vri/3)] 2 [w — exp(— 2-7ri/3)] 

+ (w-l) 2 [g 23 [w-exp(-2iri/3)] 2 +g 31 [w-exp(2iri/3)} 2 } }. (B.19b) 
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It is now clear from these formulas, together with (B.7), (B.10) and (|B.12b|) . that - 
just as in the integrable equal-coupling-constants case — also in the case with arbitrary 
coupling constants the solutions Cn( r ) feature two kinds of square-root branch points: 
those associated with triple collisions, which occur at r = t±, see (|B.10b|) . (jB.10d|) . and 
those associated with the vanishing of the derivative F'(9) = dF(8)/d8 of F(8), namely 
those associated with the 3 values 8 = 0, 8 = 27r/3, 8 = —2ir/3 (see (B.19); of course 8 is 
defined mod(27r), see (|B.6b)) h which of course correspond, as expected, to pair collisions, 
see (B.7). It is however not possible now — in contrast to the integrable case, see (jB.15|l - 
to obtain an explicit expression for the values Tb of r at which these branch points occur: 
but we infer from (B.17) and (|B.12a|) that generally there is an infinite number of such 
branch points in the complex r-plane (again, in contrast to the integrable case, when there 
only are three, see (JB.15JI ). It is on the other hand easy to obtain from these formulas an 
explicit expression for the behavior of 8(t) for r ~ r^, for instance for the branch point 
with 8{r b ) = 0, 

8{t) = (3[{t - n)/A}^ 2 + 0(\t - r b \), (B.20a) 
= i{b 12 /2) 1 ' 4 [D 2 + (r b + B)] - 1/2 . (B.20b) 
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